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I. INTRODUCTION 



The hydrophobic interaction and the excluded volume effect are two main interactions 
that determine the conformation of a polymer in a dilute solution, in space dimension d < 4r-. 
In the good solvent regime, the repulsive excluded volume effect is the dominating factor 
and the mean end-to-end distance Rn of a polymer chain with N monomers asymptoti- 
cally grows as (R%) ~ N 6 ^ d+2 \ the behavior of a self-avoiding random walk. On the other 
hand, the poor solvent regime is defined by the property that the attractive hydrophobic 
interaction between monomers dominates, where the scaling behavior is (R%) ~ N 2 l d . The 
situation is usually described by the statement that the polymer adopts a swollen confor- 
mation in a good solvent and the collapsed one in a poor solvent. The solvent where the 
repulsive and attractive interactions cancel each other is called the theta solvent, with the 
corresponding temperature being called the theta temperature, or Flory temperature, T#2. 
Tg is the temperature where the condition of the solvent changes from good to poor or vice 
versa and the collapse transition occurs. The collapse transition has been studied and the 
critical exponents have been calculated using various theoretical and computational meth- 
ods^ - —, including lattice models^— . In particular, the self-avoiding walk on a square lattice 
has been extensively studied as a model for the polymer in two dimensions, and its collapse 
transition has been studied using exact enumeration^ 1 ^ and Monte Carlo samplings*^ 1 ^ - — . 

Alternatively, phase transitions can be studied by calculating partition function zeros. 
The study of partition function zeros was initiated by the seminal paper of Yang and Lee^, 
where the zeros in the complex fugacity plane were studied to give a new insight on the 
phase transition. Subsequently the zeros in the complex temperature plane were studied by 
Fisher- 6 . With the recent advance of computational power, the study of partition function 
zeros became one the most popular methods for studying the phase transition and criti- 
cal phenomenal, and was used to study helix-coil transition of poly- alanine^ and folding 
transition of a simple model protein^. However, it was rarely used for the study of lattice 
polymers, although some preliminary qualitative results on collapse transition were reported 
for both homopolymers2S and heteropolymers^^i. 

The power of the partition function zeros method lies in its sensitivity to the onset of 
a phase transition. When the energy takes discrete values, the partition function Z is 
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expressed as 

Z = "£n{E)e^ E (1) 

E 

with n(E) being the number of states with energy E. When Z is a function of y = with 
some interaction parameter e, such as when E values are integer multiples of e, the partition 
function can be expressed in the form 



Z(y) = A(y) H(y - Vi ) (2) 

where A(y) is a function which is analytic in the whole complex plane, and y^s are the 
complex roots of the equation Z(y) = 0, called the partition function zeros. Since Z is real, 
2/iS form conjugate pairs except for the real- valued ones. By taking log and derivatives, one 
obtains the specific heat 



C N iT)= k fi\nyf 



I V _ ( V V I + fyd\ 2 ]jiA(y) 
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(3) 



where N is the size parameter of the system such as the particle number. For a system with 
the phase transition at y = y c , the locus of the zeros close in toward the positive real axis 
to intersect it at N = oo, and the singularity of Cn{T) appears in this limit. It is clear 
from Eq. (E]) that the leading behavior of such a singularity is due to the pair of partition 
function zeros closest to the real axis, called the first zeros. Therefore, by calculating the 
partition function zeros and examining the behavior of the first zeros as N — > oo, the critical 
behavior can be much more accurately analyzed than examining the behavior of Cn(T) for 
real values of the temperature, which is plagued by the noise due to the subleading terms 
containing zeros other than the first ones. 

In this work, we calculate the exact partition function zeros of the polymers on the square 
lattice up to the length N = 36, and make extrapolation of the first zero positions to estimate 
the collapse transition temperature Tq and the crossover exponent <j) (See the next section 
for the definition). The fact that our calculation is exact, along with the sensitivity of the 
partition function zeros method, allows us to estimate these quantities with reasonably high 
accuracy. 



II. THE SCALING BEHAVIOR AND THE CRITICAL EXPONENT 



The collapse transition is described by the scaling behavior of Rn near the critical tem- 
perature^, 

(R%) ~ N 2v f(rN*), (4) 
where r = (T — T e )/T e is the reduced temperature, and f(x) is a function with the property 

/(0) = 1 

\ (5) 
(x — >■ — oo), 

with /i-t being exponents that reproduce the scaling behavior of i?^ in the good and poor 
solvent regime, 

6/(d + 2) -2i/ 
^+ = 1 

2/rf-2i/ 

= — ^— ■ (6) 

In most of the studies on lattice models, the transition temperature and the critical exponents 
were usually obtained by examining the behavior of (R%) as a function of iV and T and 
fitting to Eq. (BJ, but they could also be obtained from the scaling behavior of the specific 
heaOLM: 

C N {T) ~ No+girN*), (7) 

with 



g{x) = < 



A + x~ a (x -> oo) 

const (x = 0) (8) 
A~x~ a (x -)> -oo), 



The crossover exponent <fi measures how rapidly the system undergoes the transition as the 
temperature approaches the critical temperature Tg. As will be shown later, it is directly 
related to the exponent that measures how rapidly the first zeros approach the positive real 
axis as iV — > oo. 
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III. THE MODEL 



A conformation of a polymer chain with N monomers is modeled as a two-dimensional 
self-avoiding chain of length N on a square lattice. The position of the monomer i is given 
by T{ = (k,l), where integers k and I are the Cartesian coordinates relative to an arbitrary 
origin. Chain connectivity requires |r, — rj+i| = 1, i.e., bond length is unity. Due to the 
excluded volume, there can be no more than one monomer on each lattice site, ^ Tj 
for % j. The attractive hydrophobic interaction is incorporated by assigning the energy 
— e < for each non-bonded contact between monomers. The resulting Hamiltonian is: 

^ = -e^A(r l ,r J ), (9) 



i<j 



where 



, 1 (\i — j\ > 1 and \r; — r,-| = 1) 
A(r i)r ,)=<; ' l 1 (10) 

(otherwise). 

Since the energy of the system is E = —eK where K is the number of monomer-monomer 
contacts, the partition function is expressed as a polynomial: 

A' max (A0 

Z= Yl V N (K)y K , (11) 

K=0 

where y = exp(/3e), the number of polymer conformations with contact number K, 

and K max (N) is the maximum number of possible contacts, when polymer length is N-. 

{N — 2m for m 2 < N < m(m + 1), 

" 1 J ' (12) 
N-2m-l for m(m + 1) < N < (m + l) 2 , 

where m is a positive integer. Therefore the partition function zeros can be obtained by 
enumerating the number of conformations Vl^{K) for each contact number K. The speed of 
enumeration can be increased by calculating the reduced number of conformations u^(K), 
where conformations related by rigid rotations, reflections, and translations are regarded as 
equivalent and counted only once. However, it is assumed that there is an intrinsic direction 
in the chain, so the conformations related by the exchange of labels i <H- N — i + l for all (i = 
1, ■ • • , N) are considered distinct. We note that since the rigid rotations and reflections in 
two dimensions form an eight-fold symmetry, the total number of conformations generated by 
rotations and reflections from a given two-dimensional conformation is eight. An exception 



is the straight chain, a one- dimensional conformation invariant with respect to reflection 
perpendicular to the chain. Consequently, the total number of conformations generated by 
rotations and reflections is four in this case. Therefore, the number of conformations with 
rigid rotations and reflections considered distinct, denoted by Q^(K), can be easily obtained 
by 

n N (0) = 8u N (0) -4 

n N (K) = 8u N (K) (K>Q). (13) 

Using a parallel algorithm that classifies each conformation according to the size of box it 
spans 33 , we could calculate u^{K) up to N = 36. The same quantities were calculated up 
to iV = 28 in earlier works 30,32 , which agree with the current results. 

IV. PARTITION FUNCTION ZEROS 

The partition function zeros were obtained by solving the polynomial equation 

K max {N) 

Z= n N (K)y K = 0, (14) 

K=0 

using MATHEMATICA. We observe the partition function zeros fall on a simple locus, 
more or less independent of the polymer length N (Fig. [1]). 

Although there is a relatively large gap between the positive real axis and the first zeros 
by visual inspection (Fig. [1]), the first zeros approach the positive real axis (Fig. [2]), and the 
transition temperature and the crossover exponent can be calculated from their behavior 
in the N — > oo limit. However, an oscillatory behavior is observed, due to the fact that 
there are classes of conformations whose numbers depend crucially on the parity of N. For 
example, there is only one hairpin conformation when N is even, but there are two possible 
conformations for odd iV(Fig. [3]). Therefore Ns for each parity are used separately when 
N — > oo limit is taken, so that the large error due to the oscillatory behavior is eliminated. 

The crossover exponent can be obtained from examining how fast the first zeros ap- 
proach the positive real axis as N increases^. From the scaling relation Eq. (j7J), we see that 
the partition function scales as 

\nZ N (r) ~ N^-^girN*) (15) 
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and the equation for the first zero in the first quadrant 



z(n) = o (16) 

is invariant with changing N only if 

n ~ N-+, (17) 

which is related to the corresponding complex temperature t\ as 

n = ^ (18) 

In terms of £i, Eq. ( IT71) is reexpressed as: 

ti ~ T e + const ■ (19) 
which is asymptotically equivalent to 

Vi~y c + const • (20) 

in the large N limit, where y\ = e 6 /* 1 and y c = e e//T ". From the imaginary part of Eq. (1201 

Im[ yi (iV)] ~ (21) 

the finite-size approximation of the crossover exponent is obtained: 

Mm ln{Im[y 1 (iV + 2)]/Im[ Z/l (iV)]} 

0(7Vj " ln{(iV + 2)/iV} • {ZZ) 

The expression Eq. ( 122|) reduces to the exact value of in the N — > oo limit, which is 
estimated by using the Bulirsch-Stoer (BST) extrapolation^ - — . For given m data points 
corresponding to distinct values of N, the BST extrapolation is performed by constructing 
a rational function of (l/iV) w that passes through all of these points, under the assumption 
that the leading finite-size correction is of order 0((l/iV) w ). Then, the extrapolated value 
is obtained by evaluating the function at 1/N = 0. The estimated error is defined as2&~— 

2|0_!-0_ m | (23) 

where 0_j is the value of at 1/A^ = extrapolated from the data with the 2-th point 
eliminated. The estimated error is the measure for the robustness of the extrapolated value 
with respect to perturbations in the data points, but it has no statistically rigorous confidence 



level associated with it. The estimated error can be further reduced by removing unreliable 
data obtained from N < 20, and the final result is 



(p = 0.422(12), 



(24) 



obtained from the data for even iV with 20 < N < 36. In the absence of additional 
information, we assumed the leading finite size correction to <p is of order 0(A^ _1 ) when 
performing the BST procedure, but extrapolated value of <fi does not seem to depend much 
on this assumption(data not shown). 

Once the value of <fi is determined, the transition temperature Tg can be obtained by 
estimating the point on the positive real axis where the first zeros approach in the limit of 
iV — > oo, applying the BST extrapolation procedure to the real part of Eq. (I2"0"j) . 



The resulting value of y c is 2.16(18), equivalent to Tg/e = 1.30(17), where again, the data 
for even N with 20 < iV < 36 were used. 

The finite value approximations of and y c are displayed in Figs. H] and |5] as functions 
of 1/N, along with their extrapolated values at 1/N = 0. The extrapolated value of y c in 
Fig. [5] is larger than obtained by drawing a straight line through the data points, because we 
assumed the leading behavior of y — y c being proportional to (l/iV) a422 . There is no change 
of the extrapolated value of Tg under the current precision when we use the conjectured 
exact value of the crossover exponent <fi = 3/7 14 instead of <p = 0.422. The values obtained 
in the current study are compared with those from the earlier works in Table [J Since Tg/e is 
not a universal quantity, it is displayed only for the square lattice polymer with nearest 
neighbor interaction. The maximum sizes of the polymer studied, iV max , are displayed 
wherever applicable. The results of the current study are given in the first line of the 
Table [H Although there are variations in the results reported earlier, we find that many 
of them are consistent with ours. Those that agree with our results within the estimated 
errors are indicated by boldface letters. In particular, it should be noted that the value of <j) 
obtained in the current work agrees quite well with the exact value 3/7 obtained by analytic 
calculation on the polymers on the hexagonal lattice^ 4 , which is believed to be in the same 
universality class as those on the square lattice 1 ^^* 2 ^. 



Re[ yi (N)]-y, 



c 



(25) 



S 



V. DISCUSSION 



We studied the zeros of the exact partition function of lattice polymers on square lattices 
up to chain length 36 by exhaustively enumerating the number of all possible conformations. 
We observed that the first zeros tend to approach the positive real axis as the chain length 
increases, and estimated the critical temperature Tq and the crossover exponent </> by the 
BST extrapolation. 

In contrast to Monte Carlo approaches where the calculation can be done for polymer 
lengths up to several hundreds or thousands, the chain length studied in the current study 
is much shorter, but the exactness of our data allows us to use powerful extrapolation 
methods, leading to a reasonably accurate estimation of the transition temperature and the 
crossover exponent. Furthermore, by studying the complex zeros of the partition function 
zeros, instead of examining the scaling behavior of real-valued quantities such as radius 
of gyration or specific heat, much more accurate analysis of the phase transition could be 
performed. 

It is of immediate interest to perform the exact enumeration of polymer conformations 
up to sizes where the approach of the first zeros toward the positive real axis is more visible. 
An exact enumeration has been performed using a transfer matrix for length up to 72 at 
infinite temperature^, and it would be interesting to see whether it can be generalized to 
count the number of conformations for each energy without introducing too much extra 
computational costs, in order to calculate the partition function zeros. One could also 
combine Monte Carlo methods with the partition function zeros to increase the polymer 
size, at the cost of introducing sampling error. There are indications that the locations of 
the first zeros are robust with respect to the sampling errors, a point that needs further 
investigation^. 

As a final remark, the partition function zeros method may be applied to study the 
transition behavior of heteropolymers^^i, related to the very important and interesting 
topic of protein folding. In contrast to homopolymers, the definition of large N limit is not 
so clear for a heteropolymer, so the finite size scaling argument such as the one used in the 
current study cannot be applied directly. Various methods to extract information relevant 
to the collapse and the folding transition, from the complex partition function zeros, will 
have to be explored. 
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TABLE I. The critical temperature Tg and the crossover exponent <j> obtained in the current work, 
displayed in the first line, are compared with those in the literature. Tq is displayed only for 
the model of the current work. The results that agree with ours within the estimated errors are 
indicated by boldface letters. 



Method 


lattice 




T e /e 




Exact partition function zeros 


square 


36 


1.30(17) 


0.422(12) 


Field theory 7 " 


N/A 


N/A 




X (« 0.64) 

11 v / 


Renormalization group^ 


N/A 


N/A 




£g (~ 0.86) 


Monte Carlo2 


square 


160 


1.31(6) 




Monte Carle 10 


square 


200 


1.55(15) 


0.6(1) 


Transfer matrix 11 ' 12 


square 


N/A 


1.42(4) 


0.48(7) 

V / 


Series expansion 13 


triangular 


16 




0.64(5) 


Coulomb gas method — 


hexagonal 


N/A 




| (« 0.43) 


Monte Carlo and renormalization group — 


square 


40 


1.54(7) 


0.52(7) 


Monte Carlo±£ 


hexagonal 


300 




0.5(1) 


Scanning simulation 17 


square 


240 


1.52(1) 


0.530(4) 


Recursive enrichment method 18 


square 


2048 


1.504(5) 


0.435(6) 


The pruned-enriched Rosenbluth method 19 


square 


256 


1.4993(23) 




Interacting growth walk^ - 


square 


2000 




0.419(3) 


Monte Carlo^ 


square 


1600 


1.50 


0.545(4) 


Monte Carle 22 


square 


300 


1.505(18) 




Monte Carlo 23 


square^ 


20 




0.436(7) 



A model with explicit solvent molecules. Different from the model studied in this work. 
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FIG. 1. Positions of the partition function zeros in the complex temperature (y = e^ e ) plane for 
N = 30(circles), 32(squares), 34(diamonds), and 36 (triangles). The first zeros are the ones closet 
to the positive real axis. 



13 



2.0 



1.5 



1.0 



0.5- 



0. 



ft 







Cb 



9. 




0.5 



1.0 

Re[y] 



1.5 



2.0 



FIG. 2. Positions of the first zeros in the first quadrant of the complex temperature (y 
plane for even lengths N = 10, 12, 14, • • • , 36 (circles) and for odd lengths N = 11, 13, 15, 
(squares) from left to right. The first zeros approach the positive real axis as ./V increases. 
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(a) 



(b) 
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FIG. 3. The hairpin as an example of the class of conformations whose number depends crucially 
on the parity of N. There is only one conformation for even N (a), whereas there are two possible 
conformations for odd N (b). Note that there is an instrinsic direction in a chain, indicated by an 
arrow. The dashed lines indicate the inter-monomer contacts. 
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FIG. 4. The finite size approximations of the crossover exponent, (j)(N), are shown as a function 
of 1/N for even N with TV > 20 (open circles), and the value of 4> at infinite size obtained by the 
BST extrapolation is indicated by a solid circle with an error bar. 
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FIG. 5. Values of the real part of the first zeros are shown as a function of l/N for even N with 
N > 20 (open circles), and the value for N — ^ oo obtained by the BST extrapolation is indicated 
by a solid circle with an error bar. 
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